Introdução
Aluna: Danielly Santos Severino (danyss@ufmg.br)
Orientador: Dani Gamerman (danig@ufmg.br)
Co-Orientadora: Izabel Nolau (nolau@dme.ufrj.br)
26/06/2023
Ajustar modelos para séries temporais incorporando atraso na notificação.
Verificar qualidade de ajuste dos modelos.
Avaliar erros de ajuste e erros de previsão dos modelos.
Comparar estimativas, erros de ajuste e erros de previsão entre diferentes modelos.
\(T\): Tempo atual. \(D\): Atraso máximo relevante. \(H\): Horizonte máximo de previsão.
\(n_{t,d}\): Número de eventos ocorridos no tempo \(t\) registrados após \(d\) unidades de tempo.
\(N_{t} = \sum_{d=0}^{D} n_{t,d}\): Número total de eventos ocorridos no tempo \(t\).
Por definição, existe uma estrutura de dependência entre os dados já que
\[N_{t} = \sum_{d = 0}^{D} n_{t,d}\]
E, consequentemente,
\[n_{t,0} = N_{t} - \sum_{d = 1}^{D} n_{t,d}\] para \(t = 1, ..., T, T+1, ..., T+H\).
Portanto, não é possível especificar um modelo para \(N_{t}\) e \(n_{t,d}\) incorporando todos os atrasos ao mesmo tempo respeitando a disposição dos dados.
Como \(n_{t,0}\) são as contagens sem nenhum atraso, pode ser mais interessante especificar um modelo para \(N_{t}\) e \(n_{t,d}\) com \(d = 1,...,D\) e \(t = 1, ..., T, T+1, ..., T+H\).
Modelos hierárquicos com abordagem bayesiana:
Estrutura de atraso na notificação independente.
Estrutura de atraso na notificação conjunta.
Estimação pontual e intervalar.
Previsão:
Nowcasting: Previsão para as contagens parciais.
Forecasting: Previsão para os casos futuros.
O método de Monte Carlo via Cadeias de Markov (MCMC) foi utilizado como técnica de extração de amostras para aproximar as densidades a posteriori das quantidades de interesse. Destaca-se que no MCMC, por ser um modelo altamente não linear, precisamos contar com a ajuda de bons pontos de partida para obter estimativas e previsões adequadas.
Para obter as estimativas e previsões intervalares foi utilizado o intervalo HPD (Highest Posterior Density), uma vez que esse é o intervalo de credibilidade mais estreito que contém uma proporção específica baseado na distribuição a posteriori.
A operacionalização foi realizada via STAN e R (versão 4.3.1).
O modelo com estrutura de atraso na notificação independente (MI) foi ajustado conforme as especificações abaixo.
Modelo proposto para os dados:
\[\begin{align} & n_{t,d} \sim Poisson(\lambda_{t,d}) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; N_{t} \sim Poisson(\theta_{t}) \\ \\ & \lambda_{t,d} = \exp(\alpha_{t,d}) = \dfrac{a_{d} \; c_{d}\; f_{d}\; \exp(-c_{d}\;t)} {[b_{d} + \exp(-c_{d}\;t)]^ { f_{d} + 1} } \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \theta_{t} = \dfrac{a \;c \;f \exp(-c\;t)}{[b + \exp(-c\;t)]^ { f + 1} }\\ \end{align}\]
Para \(t = 1, ..., T, T+1, ..., T+H\) e \(d = 1,...,D\).
Onde,
\[\theta_{t} > \sum_{d = 1}^{D} \lambda_{t,d}\]
Distribuições a priori:
\[\begin{align} & a_{d} \sim Gama(0.1, 0.1) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; a \sim Gama(0.1, 0.1) \\ \\ & \exp(b_{d}) \sim Normal(0, \sqrt{20}) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \exp(b) \sim Normal(0, \sqrt{20}) \\ \\ & c_{d} \sim Gama(2, 9) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; c \sim Gama(2, 9) \\ \\ & f_{d} \sim Gama(0.01, 0.01) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; f \sim Gama(0.01, 0.01) \\ \end{align}\]
O modelo com estrutura de atraso na notificação conjunta (MC) foi ajustado conforme as especificações abaixo.
Modelo proposto para os dados:
\[\begin{align} & n_{t,d} \sim Poisson(\lambda_{t,d}) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; N_{t} \sim Poisson(\theta_{t}) \\ \\ & log(\lambda_{t,d}) = \alpha_{t} + \beta_{d} \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \theta_{t} = \dfrac{a_{\theta}\; c_{\theta}\; f_{\theta} \exp(-c_{\theta}\;t)} {[b_{\theta} + \exp(-c_{\theta}\;t)]^ { f_{\theta} + 1} } \\ \\ & \exp(\alpha_{t}) = \dfrac{a_{\alpha}\; c_{\alpha}\; f_{\alpha} \exp(-c_{\alpha}\;t)} {[b_{\alpha} + \exp(-c_{\alpha}\;t)]^ { f_{\alpha} + 1} } \\ \\ & \beta_{d} = \gamma d \\ \end{align}\]
Para \(t = 1, ..., T, T+1, ..., T+H\) e \(d = 1, ..., D\).
Onde,
\[\theta_{t} > \sum_{d = 1}^{D} \lambda_{t,d}\]
Distribuições a priori:
\[\begin{align} & a_{\alpha} \sim Gama(0.1, 0.1) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; a_{\theta} \sim Gama(0.1, 0.1) \\ \\ & \exp(b_{\alpha}) \sim Normal(0, \sqrt{20}) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \exp(b_{\theta}) \sim Normal(0, \sqrt{20}) \\ \\ & c_{\alpha} \sim Gama(2, 9) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; c_{\theta} \sim Gama(2, 9) \\ \\ & f_{\alpha} \sim Gama(0.01, 0.01) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; f_{\theta} \sim Gama(0.01, 0.01) \\ \\ & \gamma \sim Normal(0,100) \\ \end{align}\]
As distribuições a priori foram especificadas de acordo com as estimativas obtidas através de uma aplicação realizada anteriormente e apresentada em Bastos et al, 2019.
Umas das formas de comparar os modelos é contrastar os erros de estimação e previsão. Para isso, foram calculadas as métricas MAE (Mean Absolute Error) e RMSE (Root Mean Squared Error).
MAE:
\[\frac{\sum_{i = 1}^{n} |\hat{y_{i}} - y_{i}|}{n} \]
RMSE:
\[\sqrt{ \frac{\sum_{i = 1}^{n} (\hat{y_{i}} - y_{i})^2}{n} }\]
O MAE e RMSE medem a magnitude média dos erros em um conjunto de estimações ou previsões.
O RMSE dá um peso maior para erros grandes. Além disso, o RMSE sempre será maior ou igual ao MAE e quanto maior a diferença entre eles, maior a variância dos erros individuais nos dados. Se MAE = RMSE, então todos os erros são da mesma magnitude.
Erros de estimação ou previsão para \(n_{t,d}\): \(y_{i} = n_{t,d}\) e \(\hat{y_{i}} = E(\lambda_{t,d}|n_{t,d})\)
Erros de estimação ou previsão para \(N_{t}\): \(y_{i} = N_{t}\) e \(\hat{y_{i}} = E(\theta_{t}|N_{t})\).
Estrutura dos casos reais de dengue no Rio de Janeiro entre janeiro de 2011 e abril de 2012.
As notificações dos casos foram feitas semanalmente e foram avaliadas para 35 semanas no total, com no máximo 10 semanas de atraso.
As contagens em vermelho foram omitidas para ajustar os modelos e verificar a qualidade das estimativas e previsões obtidas.
Em Bastos et al, 2019 \(\beta_{d}\) depende de \(\beta_{d-1}\). Aqui, adicionamos uma restrição de que que \(\beta_{d}\) é uma função linear dos atrasos.
Bastos et al, 2019 define que a estrutura conjunta de atraso na notificação apresenta um mecanismo onde as contagens dependem dos valores observados para o atraso anterior.
Aqui, adicionamos uma restrição de que os atrasos apresentam uma estrutura linear.
Observando os gráficos, parece que realmente há uma estrutura de decaimento no log com os atrasos nas notificações. Isso reforça a suposição de linearidade entre os atrasos.
MI: Modelo com estrutura de atraso na notificação independente
MC: Modelo com estrutura conjunta de atraso na notificação
BASTOS, Leonardo S., et al. A modelling approach for correcting reporting delays in disease surveillance data. Statistics in Medicine, 2019.
CHAI, T., DRAXLER, R.R. Root mean square error (RMSE) or mean absolute error (MAE)? – Arguments against avoiding RMSE in the literature. Geosci. Model Dev., 2014.
GAMERMAN, D. e LOPES, H. F. Markov Chain Monte Carlo: stochastic simulation for Bayesian inference. CRC press, 2006.
GELMAN, A., et al. Bayesian Data Analysis. 2013.
NOLAU, Izabel, GAMERMAN, Dani, BASTOS, Leonardo S. A modeling approach to forecasting data with reporting delay. Simpósio Nacional de Probabilidade e Estatística, 2022.